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New experiments, exploring the ultra-high energy tail of the cosmic ray spectrum with unprece- 
dented detail, are exerting a severe pressure on extensive air shower modeling. Detailed fast codes 
are in need in order to extract and understand the richness of information now available. Some 
hybrid simulation codes have been proposed recently to this effect (e.g., the combination of the tra- 
ditional Monte Carlo scheme and system of cascade equations or pre-simulated air showers) . In this 
context, we explore the potential of SENECA, an efficient hybrid tridimensional simulation code, 
as a valid practical alternative to full Monte Carlo simulations of extensive air showers generated 
by ultra-high energy cosmic rays. We extensively compare hybrid method with the traditional, but 
time consuming, full Monte Carlo code CORSIKA which is the de facto standard in the field. The 
hybrid scheme of the SENECA code is based on the simulation of each particle with the tradi- 
tional Monte Carlo method at two steps of the shower development: the first step predicts the large 
fluctuations in the very first particle interactions at high energies while the second step provides 
a well detailed lateral distribution simulation of the final stages of the air shower. Both Monte 
Carlo simulation steps are connected by a cascade equation system which reproduces correctly the 
hadronic and electromagnetic longitudinal profile. We study the influence of this approach on the 
main longitudinal characteristics of proton-induced air showers and compare the predictions of the 
well known CORSIKA code using the QGSJET hadronic interaction model. 



PACS numbers: 96.40.Pq,96.40.-z,13.85.-t 
I. INTRODUCTION 

Since the very first observations, ultra-high energy cos- 
mic rays (UHECR) have been an open question and a 
priority in astroparticle physics. Their origin, nature 
and possible acceleration mechanisms are still a mys- 
tery. In the last decades many experiments such as Vol- 
cano Ranch [113 Haverah Park S Yakutsk [J, Fly's 
Eye HiRes and AGASA have contributed 

for the study of UHECR's, setting up the existence of 
such high energy particles. Shorty, the Pierre Auger Ob- 
servatory will begin to explore them in unprecedented 
detail 

Due to the very low flux of high energy cosmic rays, 
measuring extensive air showers (EAS) is the only pos- 
sible technique to learn about the shape of the UHECR 
spectrum and their chemical composition. Two differ- 
ent ways have been historically applied to observe and 
analyze EAS's: ground array of detectors and optical 
detectors. Surface detectors measure a lateral density 
sample and trigger in coincidence when charged parti- 
cles pass through them. Optical detectors (i.e., fluo- 
rescence detectors) observe the longitudinal proflle evo- 
lution by measuring the fluorescence light from atmo- 
spheric nitrogen excitation produced by the ionization of 
the secondary charged particles (essentially electrons and 
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positrons). The combination of shower observables (such 
as lateral density, the depth of maximum shower devel- 
opment (Xniax) and number of muons (iV^) at detector 
observation level) and simulation techniques is the cur- 
rent way to obtain information about the primary energy, 
composition and arrival direction. For this purpose the 
shower simulation should provide all possible, and ideally 
the necessary, information to interpret measurements of 
shower parameters. We suggest the reference which 
is an interesting summary of experimental results from 
highest energy cosmic ray measurements, focused on data 
and analyzes that became available after 1999. 

Many modern shower simulation packages have been 
developed over the years. Most of them are based on the 
Monte Carlo method and simulate complete high energy 
showers with well described fluctuations in the flrst par- 
ticle interactions and realistic distributions of energy of 
shower particles. Unfortunately, the calculation of the 
gigantic showers induced by cosmic rays with energies 
above 10^* eV becomes a very difficult technical prob- 
lem. This is due to the huge number of created sec- 
ondary particles that have to be followed in the Monte 
Carlo method. As a consequence the direct simulation 
of the shower following each individual particle becomes 
practically impossible and the computation time takes 
too long. 

Recently, different ways of calculating the air shower 
development have been proposed 13, 14J. Most of them 
combine the traditional Monte Carlo scheme with a sys- 
tem of electromagnetic and hadronic cascade equations. 
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In a new one dimensional approach [15| pre-simulated 
pion-induced showers are described with parametriza- 
tions and then are superimposed to pion and kaon par- 
ticles after their first interaction points arc simulated by 
the Monte Carlo method. 

In the present work we analyze extensively the results 
obtained by the SENECA [3 code. The SENECA simu- 
lation approach is based on the Monte Carlo calculation 
of the first and final stages of the air shower develop- 
ment, and on a cascade equation system that connects 
both stages reproducing the longitudinal shower devel- 
opment. We explore mainly the fast air shower gen- 
eration for different primary energies. As an applica- 
tion of this approach we investigate the main longitu- 
dinal shower characteristics of proton, iron nucleus and 
gamma initiated air shower s up to ultra-high energy, as 
predicted by the QGS JET 0^ ^ hadronic interaction 
model. We compare SENECA results with the well tested 
CORSIKA (COsmic Ray Simulations for KAscade) sim- 
ulation code 0. 

This paper is structured as follows. In Sec^Jwe briefly 
describe the hybrid approach and the air shower mod- 
elling. In Sec. mil the method is applied to study the 
longitudinal development and the shower observables as 
well. Our main intention is to compare our predictions to 
calculations performed with the CORSIKA shower gen- 
erator, in order to verify the stability of the hybrid code 
and the reliability of the its results, in the sense of be- 
ing useful to several experimental applications. We have 
made several calculations for gamma, iron nucleus and 
proton initiated showers with different energy thresholds, 
using the thinning procedure for the Monte Carlo simula- 
tion scheme. Average and fluctuation values of ATmax and 
•Smax, their correlations and distributions are presented. 
The number of electrons (Ne) and muons {Nf^) have been 
analyzed at different observation levels. Also, we show a 
comparison between CORSIKA and SENECA CPU time 
requirements for proton showers at different primary en- 
ergies. Section HH summarizes our results. 



II. SENECA APPROACH AND THE AIR 
SHOWER MODELLING 

The main goal of this approach is the generation of 
EAS's in a fast manner, obtaining the correct descrip- 
tion of the fluctuations in showers and giving the aver- 
age values for the shower characteristics. The detailed 
technical information about the SENECA package can 
be obtained in JL6., 20J . Even though the SENECA code 
describes both longitudinal and lateral air shower devel- 
opments, the simulation scheme is used here to gener- 
ate large statistics of longitudinal shower profiles appli- 
cable mainly to the present fluorescence detectors, such 
as Pierre Auger Observatory and HiRes 0, 112 , as 
wefl to the future telescope EUSO pj . 

For the present work we track explicitly every particle 
with energy above the fraction /=i?o/l-000, where Eq is 



the primary shower energy, studying in detail the initial 
part of the shower. All secondary particles with energy 
below the mentioned fraction are taken as initial condi- 
tions to initialize a system of hadronic and electromag- 
netic cascade equations (see the suggested references). 
We use the cascade equations up to several minimum 
electromagnetic energy thresholds of 1, 3.16, 10, 31.6 and 
100 GeV for the electromagnetic component (E^f^) and 
-^min=10^ GeV for the hadronic component. Hadronic 
and electromagnetic particles with energies below the 
E^^ and E^^ are no longer treated by the cascade 
equations but traced again in the Monte Carlo scheme. 
The hadronic interactions at high energies are calculated 
with the QGSJETOl ^3 model while the interactions at 
low energies with GHEISHA 24]. The electromagnetic 
interactions are treated by the EGS4 code. The 
adopted kinetic energy cutoffs for all simulations were 
50 MeV (0.3 MeV) for hadrons and muons (electrons 
and positrons). The forementioned energy thresholds 
were chosen in order to verify qualitatively the depen- 
dence and performance of the cascade results, trying to 
minimize the CPU time without the introduction of sig- 
nificant errors. All SENECA simulations were performed 
with 1.2.2. version. 



III. RESULTS AND COMPARISONS 

In this session we apply the SENECA simulation ap- 
proach to generate gamma, iron nucleus and proton in- 
duced air showers at fixed primary energies and explore 
the longitudinal development. Although the simulation 
of showers at fixed energies is not a very realistic applica- 
tion we intend in the present work to compare quantita- 
tively SENECA and CORSIKA results. One important 
reason for this comparison is to optimize the compro- 
mise between simulation time usage and accuracy in the 
description of fluctuations. In this particular case, such 
detailed study can be useful to many experiments which 
use the fluorescence technique. 

In Fig. ^ we illustrate individual longitudinal proflles 
generated by SENECA (hybrid approach) code. Dis- 
played cases correspond to median longitudinal profile 
produced for 1,000 gamma (dashed), iron nucleus (dot- 
dashed) and proton (solid) induced showers at 10"^^ eV, at 
zenith angle 6=45° and free first interaction point. The 
shaded bands which follow each median line correspond 
to the upper limit of 68% of confidence level. 

In order to make a simple comparison Fig.[21illustrates 
the upper limit of 68% of confindence level for 1,000 
SENECA profiles with 10 random showers simulated with 
CORSIKA. It is possible to see a very reasonable agree- 
ment among the several longitudinal profiles. In spite 
of that, we can verify that CORSIKA produces, appar- 
ently, more fluctuations related to individual proton and 
gamma air showers developments. All simulations per- 
formed by the CORSIKA code were obtained by using 
the thinning factor if=10~® and the same energy thresh- 
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FIG. 1: Median longitudinal profile of electrons-positrons for 
10 random gamma (dashed), iron nucleus (dot dashed) and 
proton (solid) induced showers at primary energy of 10^^ eV, 
with zenith angle 0=45°, calculated by the SENECA scheme, 
by using the QGSJETOl hadronic interaction model. The 
shaded bands which follow each median line represent the 
upper limit of 68% of confidence level, for 1,000 simulated 
events. 



old of SENECA simulations, 0.3 MeV (50 MeV) for pho- 
tons and electrons (hadrons and muons). 

To analyze more carefully the longitudinal profiles pro- 
duced with the hybrid approach, we illustrate in Fig. O 
the number of electrons-positrons generated by 1,000 pro- 
ton initiated showers at energy of 10^^ eV, with zenith 
angle 45°, at arbitrary observation depths of (a) 200, (b) 
400, (c) 600, (d) 800, (e) 1,000, and (f) 1,200 g/cm^. We 
compare the hybrid results (solid line) to the results ob- 
tained with CORSIKA (dashed line) . The expected fluc- 
tuations due to very first proton interactions are reflected 
on panel (a) and (b) , producing a difference of ~ 4% and 
~ 3% between the distribution mean values. Panels (c) 
and (d) show discrepancies for the mean values of about 
7% and 6%. In contrast to this, the discrepancies basi- 
cally disapear at panels (e) and (f), both obtaining an 
agreement of about 99.5%. In all found discrepancies, 
SENECA produces more particles than CORSIKA. Ac- 
cording to some discrepancies are expected in the 
values produced for some shower quantities by the hy- 
brid scheme and, apparently, are due to the combination 
of the used binning of discrete energy in the numerical 
solutions and the minimum energy threshold used for the 
electromagnetic cascade equations. 

Such differences between the results of both codes at 
the first stages of the EAS development should not be 
crucial since depths around and after the shower max- 
imum are the most important ones for all fluorescence 
experiments. 
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FIG. 2: Longitudinal profile of electrons-positrons for 10 ran- 
dom gamma, iron nucleus and proton initiated air showers at 
primary energy of 10^^ eV, with incident zenith angle of 45°, 
calculated by the CORSIKA (long-dashed) code, with the 
QGSJETOl hadronic interaction model. The shaded curve 
in each panel represent the upper limit of 68% of confindence 
level for 1,000 events simulated with the SENECA code. 



A. Analyses of influence of input parameters 
on shower quantities 

In order to verify the dependence of the SENECA re- 
sults, as suggested by the SENECA authors, on the com- 
bination of used energy binning and minimum energy 
threshold (E^^) we have simulated air showers induced 
by gamma and proton primaries, with different initial 
conditions. Such primary particles have been chosen due 
to the different shower development in the atmosphere. 
The simulation development for gamma-induced showers 
deals mostly with the electromagnetic processes while the 
simulation development for proton (or any other hadron) 
showers deals with hadronic and electromagnetic pro- 
cesses. In other words, we are going to have different 
request of the cascade equation system for proton and 
gamma primary particles in the hybrid simulation. 



Parameter 



For the test with protons, 500 showers were generated 
of primary energy lO""^^ eV, with incident zenith angle of 
45°, calculated with the hybrid method using 30, 40 and 
50 bins (10 bins) in the numerical solutions of electro- 
magnetic (hadronic) cascades, and different minimum en- 
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FIG. 4: Distribution of S'max normalized by the primary en- 
ergy in GeV. Results are shown for 500 proton showers of 
energy 10^^ eV, with incident zenith angle of 45°, calcu- 
lated with the hybrid scheme with different minimum energy 
thresholds and binnings of discrete energy, for the electro- 
magnetic cascade equations. The thick solid line refers to the 
CORSIKA distribution. Both simulation codes generated the 
results by using the QGSJETOl hadronic interaction model. 
Panels (a), (b) and (c) refer to the parameter values of 30, 40 
and 50 bins of discrete energy, respectively, in the numerical 
solutions of the hybrid simulations. 
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FIG. 3: Shower distribution in number of electrons at dif- 
ferent slanth depths. Results are shown for 1,000 showers, 
at 45°, generated by primary protons of energies 10^^ eV 
calculated with the SENECA method (solid) and CORSIKA 
code (dashed), both using the QGSJETOl hadronic interac- 
tion model. Each panel represents a particular and arbitrary 
slanth depth: (a) 200, (b) 400, (c) 600, (d) 800, (e) 1,000 and 
(f) 1,200, all in g/cm^ 



ergy thresholds for the electromagnetic cascade equations 
E^f^= 1, 10 and 100 GeV. Following we adopted 
£:had^ 104 GeV. 

Fig. 21 shows the distribution of S'max normalized by 
the primary energy in GeV, calculated with 500 primary 
protons by using 30 (panel a), 40 (panel b) and 50 (panel 
c) bins in the numerical solutions. The solid, dashed and 
dot-dashed lines correspond to E^^= 1, 10 and 100 GeV, 
respectively, while the thick solid line refers to the distri- 
bution obtained with the CORSIKA code. 



Although we have modified the inputs for the SENECA 
simulations, the average values of Smax/E do not vary 
significantly and are in total agreement (99% on the av- 
erage) with the average Smax/E obtained with the COR- 
SIKA shower generator. However, it is possible to verify 
in Fig. ^ that the S'max distributions for the hybrid code 
are visibly wider for the minimum energy threshold val- 
ues of 10 and 100 GeV, when compared to CORSIKA, 
showing a different behaviour for the Smax fluctuations in 
both codes. Besides, such fluctuation description changes 
dramatically the characteristic asymmetrical shape of the 
Smax distribution. The average (fluctuation) values of 
SraaxI E Calculated with the hybrid scheme are presented 
in Table and correspond to the distributions illustrated 
in Fig.^ The average value obtained with the CORSIKA 
code is 0.646 while the width distribution is 2.06x10^'^. 

Such Smax fluctuation descriptions for the mini- 
mum electromagnetic energy threshold values of 10 and 
100 GeV can clearly be seen as very strong dependence 
of shower quantities on initial parameters and are much 
larger when compared to the standard deviation value 
predicted by the CORSIKA code. One equally impor- 
tant aspect is that we have verified the dependence of 
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TABLE L Average values of Smax (standard deviation), normalized by the primary energy in GeV, obtained by hybrid sim- 
ulations of 500 proton showers of primary energy E = 10^^ eV, with incident zenith angle of 45°. The predictions refer to 
the hybrid scheme by using different minimum electromagnetic energy threshold and distint binnings of discrete energy. The 
average value obtained with the CORSIKA code is 0.646 while the sigma is 2.06x10"^. 
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FIG. 5: Distribution of 5max normalized by the primary en- 
ergy in GeV. Results are shown for 500 primary gamma show- 
ers of energy 10^® eV, with zenith angle of 45° calculated 
with the SENECA and CORSIKA (thick solid hne) schemes 
by using the QGSJETOl hadronic interaction model. The 
SENECA simulations were performed with different minimum 
energy thresholds for the electromagnetic cascade equations 
and the input value of 30 bins in the numerical solutions. 



FIG. 6: Median values of the shower size at maximum air 
shower development (Smax), normalized by the primary en- 
ergy in GeV, related to proton-induced showers as a function 
of primary energy. The solid (long-dashed) line represents a 
band containing 68% of confindence level for 500 (100) events 
generated by the SENECA (CORSIKA) code, at S = 45°, 
using the QGSJETOl hadronic interaction model. The full 
square (full circle) represents the median values of Smax ob- 
tained with the SENECA (CORSIKA) code. 



quantities on higher values of Ef^^, while 0| expects 
the dependences on lower values. 

The same unsual fluctuation description for S'max can 
also be visualized in Fig.El which shows the 5max/E dis- 
tribution calculated with primary gammas by using 30 
bins as input value in the numerical solutions. The thick 
solid line refers to the distribution obtained with the 
CORSIKA code. The average (sigma) values obtained 
from SENECA are 0.7 (4.xl0-2), 0.698 (4.02 xlO^^) and 
0.694 (7.79x10-2), for E^^= 1, 10 and 100 GeV, re- 
spectively. The CORSIKA predicts the average S'max/E 
distribution value of 0.694 and the width distribution of 
4.74x10-2. 

Hence, these checks basically verify that, apparently. 



the minimum electromagnetic energy thresholds of 10 
and 100 GeV may produce artificial fluctuations on S'max- 
Moreover, such fluctuation description can introduce sta- 
tistical errors on analyses of event-to-event and could be 
possible critics in pratical applications of SENECA in 
fluorescence event reconstruction. 

The increase in fluctuations with increasing i?™;;, is 
difhcult to understand since the larger the E^^, the 
greater the weighted of the Monte Carlo portion of the 
hybrid scheme. 

Fig. El shows the median shower size values, normal- 
ized by the primary energy in GeV, at the depth of 
shower maximum. Each band contains 68% of confidence 
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FIG. 7: Distribution of the depth of maximum air shower de- 
velopment shown for 1,000 proton, iron nucleus and gamma 
showers at a particular primary energy of 10^^ eV, with inci- 
dent zenith angle of 45° , calculated with the hybrid technique. 



level and is related to proton initiated showers as a func- 
tion of primary energy (i?>10^^ eV). The solid (dashed) 
line illustrates the predictions of the SENECA technique 
(CORSIKA), using the QGSJETOl model. We have used 
^m?n= 1 GeV and 30 bins of discrete energy as the input 
SENECA parameters. 

The values predicted by both codes have a visible 
agreement (> 99%). SENECA tends to produce sys- 
tematically higher values of S'max than CORSIKA, up to 
-|^q19.5 g^j-^j lower values for energies >10^*' eV. Also, it 
seems that SENECA produces more fluctuations related 
to smaller shower sizes at energies around 10^*^ eV. 



2. X„ 



Parameter 



As discussed in Sec. U fluorescence detectors observe 
the longitudinal profile evolution by measuring the fluo- 
rescence light. The greatest advantage of this technique 
is that it allows the estimation of the number of charged 
particles as a function of depth in the atmosphere, from 
the measured data [2^ |23, IH, , which makes possi- 
ble, in turn, the estimation of the depth of maximum 
development of the shower, Xmax- 

In principle, obtaining the values of ATmax and/or 5max, 
by the fluorescence technique, and their respective fluctu- 
ations, by Monte Carlo, one should be able to reconstruct 
the shower energy and infer the identity of the primary 
cosmic ray |^ |22| . Fig[7| illustrates the potential of the 
Amax distribution, generated with the hybrid scheme, to 
distinguish possible primary signatures, as calculated us- 
ing SENECA. 

We have also checked the behaviour of average X^ax 
values for gamma, iron nucleus and proton induced air 
shower. 



Fig. |51 shows the average depth (bottom panel) and 
fluctuation (top panel) values of shower maximum 
(Amax) for the energy range of 10^^-10^°-^ eV, with in- 
cident zenith angle oi 9 — 45°. All simulations per- 
formed in this plot were generated with the QGSJETOl 
hadronic interaction model. The long-dashed line rep- 
resents 500 showers induced by gamma (diamond), iron 
nucleus (square) and proton (circle) primaries at each 
energy generated by SENECA. The CORSIKA predic- 
tions correspond to 100 simulated showers and are repre- 
sented by solid lines and full symbols 30] . Only in Fig.|H| 
the CORSIKA results obtained for gamma showers were 
generated by using the PRESHOWER code. All other 
CORSIKA simulations in the present work were calcu- 
lated without the PRESHOWER code. 

The values predicted by both codes are in very good 
agreement, showing almost identical values for A^ax and 
fluctuations, for all primary particles. One interesting 
aspect to be mentioned is the strong influence of the ge- 
omagnetic fleld on depth of maximum development in 
gamma showers, generated by CORSIKA code |39| im- 
plemented with the PRESHOWER program [sj , at ener- 
gies above 10^" eV. The gomagnetic field decelerates the 
air gamma shower development, modifying the depth in 
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FIG. 8: Average depth (bottom panel) and sigma (top panel) 
values of maximum air shower development (Xmax) related to 
gamma, iron nucleus and proton induced showers as a func- 
tion of extremely high primary shower energy, at 6 — 45°, 
using the QGSJETOl hadronic interaction model. The long- 
dashed lines represent 500 events per energy induced by 
gamma (diamond), iron nucleus (square) and proton (circle) 
primaries generated by the SENECA. The CORSIKA predic- 
tions correspond to 100 simulated showers and are represented 
by solid lines and full symbols [SOi] . 
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FIG. 9: Distribution of the depth of maximum air shower de- 
velopment. Results are shown for 500 proton showers of en- 
ergy 10^^ eV, with 45°, calculated with the SENECA by using 
different minimum energy thresholds and binnings of discrete 
energy, for the electromagnetic cascade equations. The thick 
solid line refers to the CORSIKA distribution. Both simu- 
lation codes generated the results by using the QGSJETOl 
hadronic interaction model. Panels (a), (b) and (c) refer to 
the parameter values of 30, 40 and 50 bins of discrete energy, 
respectively, in the numerical solutions of the hybrid simula- 
tions. 



which occurs the maximum development. The average 
(fluctuation) value of X^ax decreases from 1,213 (314) 
g/cm^ at 10^0 eV to 990 (35) g/cm^ at lO^"-^ eV. Altough 
we have considered the geomagnetic field in our simula- 
tions, it seems that SENECA do not consider the inter- 
action of ultra high energy photons with the Earth's geo- 
magnetic field before entering the Earth's atmosphere. 
Such effect absence produces the discrepancy seen in 
gamma curves at energies above and around 10^" eV. Ref- 
erences (33,113 discuss very high energy gamma showers. 

Fig. 1^ confirms the strength of the agreement by 
showing the distribution of Xmax values, produced by 
500 proton showers, for the particular primary energy 
10^^ eV, with incident zenith angle of 45°, calculated 
with SENECA and CORSIKA (thick sofid fine) codes 
and QGSJETOl hadronic interaction model. We used 
the same input parameters for SENECA as in Fig. 0] 
The mean (fluctuation) Xmax values calculated with 
the hybrid scheme are presented in Table ^1 These 
values are in total agreement with CORSIKA values, 
(-''^max)=774 g/cvc? (cr=65 g/cm^), and show, apparently, 
that the influence of minimum energy threshold and bin- 
ning of discrete energy on X^ax values is negligible in 
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FIG. 10: Distribution of the depth of maximum air shower 
development. Results are shown for 500 gamma showers of 
energy 10^^ eV, with incident zenith angle of 45°, calculated 
by the hybrid scheme with different minimum energy thresh- 
olds for the electromagnetic cascade equations. The thick 
solid line refers to the CORSIKA distribution. The obtained 
values refer to 30 bins of discrete energy in the numerical 
solutions of the hybrid simulations. 



hadron initiated showers. 

We have made the same X^ax verification for gamma 
showers in Fig. ^] and we verify that CORSIKA and 
SENECA distributions are in good agreement. SENECA 
predicts the average (fluctuation) Xmax distribution val- 
ues of 981 g/cm^ (81 g/cm^), 989 g/cm^ (85 g/cm^) and 
985 g/cm^ (101 g/cm^), for minimum electromagnetic 
energy thresholds of 1, 10 and 100 GeV, respectively. 
CORSIKA code, without PRESHOWER, produces the 
average and width distribution values of 976 g/cvc? and 
97 g/cm^, respectively. Due to the extreme fluctua- 
tions in shower development induced by primary gam- 
mas, which affect directly the position of shower max- 
imum development, the average and fluctuation results 
obtained by SENECA seem to be reasonable, in spite of 
the fact that again the fluctuation values increase with 
the electromagnetic minimum energy thresholds of 10 
and 100 GeV, over those of CORSIKA. 



3. Srr 



-Xn 



Correlation 



We also analyzed the correlation between iSmax and 
Xniax, which are important longitudinal shower quanti- 
ties on event reconstruction. We compared SENECA re- 
sults with CORSIKA predictions for gamma and proton 
induced air shower. 

Fig.llllshows the correlation between these parameters 
for both simulation schemes. We simulated 500 gamma 
showers for each code, of energy 10^^ eV, incident zenith 
angle 0=45° and free flrst interaction point. Hybrid re- 
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TABLE 11: Average (standard deviation) values of Xmax, all in g/cm^, obtained by hybrid simulations of 500 proton showers 
of primary energy E = 10^^ eV, with incident zenith angle of 45°. The predictions refer to the hybrid scheme by using 
different minimum electromagnetic energy threshold and distint binnings of discrete energy. The average value obtained with 
the CORSIK A code is 774 g/cm^ while the sigma is 65 g/cm^ . 
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^max (o" ) 
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775 (63) 


776 (64) 


£,X=10 GeV 


775 (62) 


772 (61) 


774 (65) 


£;™'„=100 GeV 


774 (69) 


771 (66) 


777 (70) 



suits are shown for minimum enere;y thresholds 1 

0.7 111111 

(crosses) and 100 GeV (squares), and 30 bins of discrete 
energy in the numerical solutions. The full circles denote 
the corresponding values for CORSIKA. It is possible to 
verify the existence of large fluctuations in gamma air 
showers at this particular primary energy. 

The small box in the figure encloses, approximately, 
the highest density of correlation points for both codes: 
61% and 59% of the total number of events simulated by 
SENECA, for 1 and 100 GeV, respectively, and 

64% of the total showers generated with CORSIKA. 
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FIG. 11: The correlation between Xmax and Smax/ for gamma 
induced showers at primary energy 10^^ eV, at zenith angle 
f = 45°, obtained with SENECA and CORSIKA codes. The 
full circles represent 500 showers simulated with CORSIKA, 
while square and cross symbols illustrate 500 events generated 
with the hybrid method, by using minimum electromagnetic 
energy thresholds E^^= 1 and 100 GeV, respectively. 

The number of gamma showers predicted by SENECA 
and CORSIKA vary considerably from one quadrant to 
another in Fig. 1111 indicating qualitative and quantitative 
differences in the character of the fluctuations in both 
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FIG. 12: The correlation between Xmax and Snuix/E for pro- 
ton induced showers at primary energy 10^^ eV, at zenith 
angle 6 = 45°, obtained with SENECA and CORSIKA 
codes. The full circles represent 500 showers simulated with 
CORSIKA, while square and cross symbols correspond to 
500 showers generated with the hybrid method, by using 
minimum electromagnetic energy thresholds E'j^^= 1 and 
100 GeV, respectively. 



codes. 

The fractions of events inside regions I, II, III and IV 
(see the figure) are, respectively, 10% (10%), 24% (27%), 
3% (2%) and 2% (2%), for the minimum energy threshold 
of 1 (100) GeV. In the same regions, the corresponding 
CORSIKA fractions are 11%, 17%, 6% and 2%. 

Notoriously, CORSIKA shows in region I a par- 
ticular structure related to maximum depth develop- 
ment (Xniax>940 g/cm^) of particularly large showers 
{Smax/E >0.7G5). Such structure is absent in SENECA. 

On the other hand, the hybrid scheme predicts a sig- 
nificantly number of ordinary and small showers in region 
II (118 and 137 showers for E^^^= 1 and 100 GeV, re- 
spectively), which achieve higher values for the depth of 
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FIG. 13: Average values of Xmax and Smax/ for 1,000 proton 
induced showers at primary energy 10^^ eV, at zenith angle 
f = 45°, obtained with SENECA and CORSIKA (thick solid 
line) codes (see the text). 



FIG. 14: Fluctuation of the position of shower maximum and 
maximum shower size, for 1,000 proton induced showers at 
primary energy 10^" eV, at zenith angle = 45°, obtained 
with SENECA and CORSIKA (thick solid line) codes (see 
the text). 



maximum development, a behavior that is not as strong 
in CORSIKA (84 showers). 

The same kind of comparison is presented in Fig. 1121 
for proton showers. The same input parameters as in the 
previous figure were used to simulate the three sets of 
500 proton showers each displayed in Fig. 

It is apparent from the figure that the hybrid approach 
generates a more scattered distribution of showers, with 
wider tales, than CORSIKA does. 

This is confirmed by counting the number of events 
inside the small box shown in the figure. SENECA ex- 
pectation amounts to 80% (69%) of the proton showers 
falling inside the box for E^^— 1 and 100 GeV, while 
CORSIKA expectation is 85%. 

We can verify in SENECA predictions that the in- 
put parameter £'^jj=100 GeV produces larger fluctua- 
tions related to proton shower size {o'snia.x/E)^ describ- 
ing a sort of symetrical distribution of events outside the 
box. The standard deviation of proton shower size, ob- 
tained for this particular input parameter in SENECA, is 
3.62x10^^ while the corresponding value for CORSIKA 
is 1.89x10^^. Such symetrical statement is confirmed in 
Fig.EI 

Figs. El and illustrate the dependence of the av- 
erage values, and their fluctuation, of X^ax and S'max 
on input parameters of electromagnetic cascade equa- 
tions solved in SENECA scheme. We have generated 
1,000 proton showers, at primary energy 10^^ eV and 
zenith angle 6 = 45°, using 30, 40 and 50 bins of pos- 
sible discrete energy for each minimum energy thresh- 
old; total of 15,000 simulated events. As it has been 



noted in top panel of Fig. ^1 the average values pro- 
duced by SENECA are in perfect agreement (>99% for 
both quantities) with the corresponding results obtained 
from CORSIKA (thick sohd line). 

However, Fig.ll4lshows that the fluctuations related to 
'S'max (top panel) have large discrepancies for a minimum 
electromagnetic energy threshold of 100 GeV. In this 
case, there is more than a factor of 2 between SENECA 
and CORSIKA expected fluctuation levels (cr ~ 4.15 x 
10~^, averaged over number of bins, vs. a ~ 2.034 x 10~^ 
in the case of CORSIKA). 

On average, the discrepancies in a between SENECA 
and CORSIKA decrease systematically at progressively 
lower values of Ef^^^. 

Even if fluctuations in X^ax can be up to 8% larger 
in SENECA than in CORSIKA (for E^^=100GeY) this 
discrepancy is of small practical relevance when trans- 
lated into units of depth (i.e., in g/cm^). 



4- - Number of muons 

To ensure the consistency of other longitudinal shower 
components in SENECA, we have also compared the 
number of muons produced by both SENECA and COR- 
SIKA codes. 

Fig. 1151 shows the muon component for 1,000 proton 
initiated showers at energy of 10^^ eV and incident zenith 
angle of 45°, at observation depths of (a) 200, (b) 400, 
(c) 600, (d) 800, (e) 1,000, and (f) 1,200 g/cm^. As we 
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can seen, both distribution agree very well. The average 
number of muons are very similar ( 99% of coincidence) 
when related to the first stages of shower development, 
panels (a), (b) and (c). At the same time, discrepan- 
cies become greater with depths closer to the sea level, 
but are on average less than ^ 3% (2% at 800 g/cm^, 
3.6% at 1,000 g/cm2 and 3.3% at 1,200 g/cm^). Such 
small discrepancies can be considered negligible, because 
they are undetectable to any experiment which measures 
the correlated number of particles. Furthermore, the 
fluctuations (width distribution) corresponding to each 
panel agree well with those predicted in CORSIKA code. 
Consequently, the muon distributions seem to be well 
described in the hybrid scheme when compared to the 
CORSIKA ones. 



5. Processing time consumption 

As a final comparison, we show in Table IIIII the im- 
pressive numbers related to the rate of time-consumption 
(Tcors/Tsen) between the CORSIKA and SENECA 
codes. All simulations generated with CORSIKA were 
obtained by using the thinning factor 

In order to compare the time-consuming of the hybrid 
method, by using the thinning procedure in the Monte 
Carlo calculation, we have generated air showers with 
thinning level of 10~^ for the Monte Carlo scheme. 

In Table lllll Rate-1 refers to the hybrid simulations 
without the thinning procedure while Rate-2 considers 
the SENECA simulations by using the thinning proce- 
dure with thinning level of 10^^ in the Monte Carlo 
scheme. 



IV. DISCUSSION AND CONCLUSIONS 

In the present work we analyzed the practical potential 
of SENECA, a very fast hybrid tri-dimensional code, for 
the simulation of the longitudinal development of exten- 
sive air showers at high energies. We take as reference 
the well known and extensively tested CORSIKA code, 
which is based on a much slower, but highly reliable, full 
Monte Carlo simulation. The QGSJETOl hadronic in- 
teraction model is used throughout the paper for both 
codes. 

Although a careful analysis of many shower quantities 
discussed here is strongly model dependent, we are confi- 
dent that the present study is able to show the potential 
and limitations of SENECA, in its present version, for a 
practical application to the analysis of fluorescence data 
on ultra- high energy cosmic rays. 

The consistency of the SENECA scheme was tested 
and it proved to be very stable for energies above 10^^ 
eV. The results obtained by both codes agree well with 
negligible discrepancies for most quantities analyzed. 

Apparently, minimum electromagnetic energy thresh- 
old values of 10 and 100 GeV may produce artificial 



150 
125 
100 
75 
50 
25 




_ 1 1 1 1 1 1 1 1 1 1 1 1 1 


1 1 1 1 1 1 1 I. 


Proton Showers HTP 




\ (f)X, =1,200 g/cm' 




- _r^ 

: , , , . 1 . r — r^^V 1 , , , 


1 , P^'— -, , .: 



2.00e+07 



3.50e+07 5.00e+07 



6.50e+07 




4.00e+03 



9.60e+06 



FIG. 15: Shower distribution in number of muons at differ- 
ent slanth depths. Results are shown for 1,000 showers, at 
45°, generated by primary protons of energies 10^^ eV cal- 
culated with the SENECA method (solid) and CORSIKA 
code (dashed), both using the QGSJETOl hadronic interac- 
tion model. Each panel represents a particular and arbitrary 
slanth depth: (a) 200, (b) 400, (c) 600, (d) 800, (e) 1,000 and 
(f) 1,200, all in g/cm^ 



fiuctuations on S'max shower quantity. These could be 
responsible for the introduction of systematic errors on 
the analysis of isolated events and should be taken into 
account when applying SENECA to fluorescence recon- 
struction. 

In any case, the undisputable bounty of SENECA is 
velocity. The rates of processing time-consumption be- 
tween CORSIKA and SENECA are, to say the least, im- 
pressive over the wide range of primary shower energy 
tested here. 
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TABLE III: Rate of time-consuming between CORSIKA and SENECA codes, i.e., Rate=TcoRs/TsEN- Rate-1 refers to the 
hybrid scheme with the traditional Monte Carlo method, while Rate-2 considers the hybrid scheme using the thinning level of 
10~® in the Monte Carlo. All CPU times used for the calculation of the rate refer to a 2.1 GHz AMD Athlon processor. 



logio^Jo [eV] Rate-1 Rate-2 
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Finally, as a word of caution, special attention must be 
paid to the selection of the several input parameters of 

SENECA in practical applications, since they may affect 
the results in a non-trivial and somehow unpredictable 
way. 
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